Applying term consistency to an equality constrained interval global optimization problem

ABSTRACT

One embodiment of the present invention provides a system that solves a global optimization problem specified by a function | and a set of equality constraints q 1 (x)=0 (i=1, . . . , r), wherein | is a scalar function of a vector x=(x 1 , x 2 , x 3 , . . . x n ). During operation, the system receives a representation of the function ƒ and the set of equality constraints and stores the representation in a memory within a computer system. Next, the system and performs an interval global optimization process to compute guaranteed bounds on a globally minimum value of the function ƒ (x) subject to the set of equality constraints. Performing this interval global optimization process involves, applying term consistency to the set of equality constraints over a subbox X, and excluding portions of the subbox X that violate the set of equality constraints.

RELATED APPLICATION

[0001] The subject matter of this application is related to the subjectmatter in a co-pending non-provisional application by the same inventorsas the instant application and filed on the same day as the instantapplication entitled, “Applying Term Consistency to an InequalityConstrained Interval Global Optimization Problem,” having serial numberTO BE ASSIGNED, and filing date TO BE ASSIGNED (Attorney Docket No.SUN-P6446-SPL).

BACKGROUND

[0002] 1. Field of the Invention

[0003] The present invention relates to performing arithmetic operationson interval operands within a computer system. More specifically, thepresent invention relates to a method and an apparatus for using acomputer system to solve a global optimization problem includingequality constraints with interval arithmetic and term consistency.

[0004] 2. Related Art

[0005] Rapid advances in computing technology make it possible toperform trillions of computational operations each second. Thistremendous computational speed makes it practical to performcomputationally intensive tasks as diverse as predicting the weather andoptimizing the design of an aircraft engine. Such computational tasksare typically performed using machine-representable floating-pointnumbers to approximate values of real numbers. (For example, see theInstitute of Electrical and Electronics Engineers (IEEE) standard 754for binary floating-point numbers.)

[0006] In spite of their limitations, floating-point numbers aregenerally used to perform most computational tasks.

[0007] One limitation is that machine-representable floating-pointnumbers have a fixed-size word length, which limits their accuracy. Notethat a floating-point number is typically encoded using a 32, 64 or128-bit binary number, which means that there are only 2³², 2⁶⁴ or 2¹²⁸possible symbols that can be used to specify a floating-point number.Hence, most real number values can only be approximated with acorresponding floating-point number. This creates estimation errors thatcan be magnified through even a few computations, thereby adverselyaffecting the accuracy of a computation.

[0008] A related limitation is that floating-point numbers contain noinformation about their accuracy. Most measured data values include someamount of error that arises from the measurement process itself. Thiserror can often be quantified as an accuracy parameter, which cansubsequently be used to determine the accuracy of a computation.However, floating-point numbers are not designed to keep track ofaccuracy information, whether from input data measurement errors ormachine rounding errors. Hence, it is not possible to determine theaccuracy of a computation by merely examining the floating-point numberthat results from the computation.

[0009] Interval arithmetic has been developed to solve theabove-described problems. Interval arithmetic represents numbers asintervals specified by a first (left) endpoint and a second (right)endpoint. For example, the interval [a, b], where a<b, is a closed,bounded subset of the real numbers, R, which includes a and b as well asall real numbers between a and b. Arithmetic operations on intervaloperands (interval arithmetic) are defined so that interval resultsalways contain the entire set of possible values. The result is amathematical system for rigorously bounding numerical errors from allsources, including measurement data errors, machine rounding errors andtheir interactions. (Note that the first endpoint normally contains the“infimum”, which is the largest number that is less than or equal toeach of a given set of real numbers. Similarly, the second endpointnormally contains the “supremum”, which is the smallest number that isgreater than or equal to each of the given set of real numbers.)

[0010] One commonly performed computational operation is to performequality constrained global optimization to find a global minimum of anonlinear objective function subject to nonlinear equality constraintsof the form q₁(x)=0 (i=1, . . . , r). This can be accomplished bydeleting boxes, or parts of boxes that do not satisfy one or moreequality constraints, or by unconditionally deleting boxes that cannotcontain the global minimum |*.

[0011] To delete part or all of a subbox that fails to satisfy one ormore equality constraints, term consistency and the interval Newtonmethod can be used in combination.

[0012] To unconditionally delete part or all of a box that cannotcontain the equality constrained global minimum, the smallest upperbound |_bar so far computed on the global minimum |* can be used.

[0013] In the presence of equality constraints, it can be difficult toidentify a feasible point x that is guaranteed to satisfy all the givenequality constraints. However, it is relatively easy to prove that a boxX contains a feasible point. Given such a box, sup(|(X)) provides avalue of |_bar. As this process progresses, |_bar is updated using|_bar=mince bar, sup(|(X))), given that the box X has been proved tocontain a feasible point.

[0014] Given an |_bar value, the |_bar criterion is applied as aninequality constraint to delete parts or all of a given subbox.

[0015] If |_bar is the smallest upper bound so far computed on |, thenany point x for which |(x)>|_bar can be deleted. Similarly, any box Xcan be deleted if in|(|(X))>|_bar.

[0016] In addition to solving for boxes with feasible points and the|_bar criterion, the John conditions can be solved to delete parts orall of subboxes and to prove that a given subbox contains a feasiblepoint. (The John conditions are described in “Global Optimization UsingInterval Analysis” by Eldon R. Hansen, Marcel Dekker, Inc., 1992.)

[0017] Solving the John conditions, works best “in the small” when theobjective function | is approximately quadratic and satisfiedconstraints are approximately linear. For large boxes containingmultiple stationary points, solving the John conditions will often notsucceed in deleting much of a given box. In this case the box is splitinto two or more sub-boxes that are then processed independently. Bythis mechanism all the equality constrained global minima of a nonlinearobjective function can be found.

[0018] One problem is applying this procedure to large n-dimensionalinterval vectors (or boxes) that contain multiple local minima. In thiscase, the process of splitting in n-dimensions can lead to exponentialgrowth in the number of boxes to process.

[0019] It is well known that this problem (and even the problem ofcomputing “sharp” bounds on the range of a function of n-variables overan n-dimensional box) is an “NP-hard” problem. In general, NP-hardproblems require an exponentially increasing amount of work to solve asn, the number of independent variables, increases.

[0020] Because NP-hardness is a worst-case property and because manypractical engineering and scientific problems have relatively simplestructure, one problem is to use this simple structure of real problemsto improve the efficiency of interval equality constrained globaloptimization algorithms.

[0021] Hence, what is needed is a method and an apparatus for using thestructure of a nonlinear objective function to improve the efficiency ofinterval equality constrained global optimization software. To this end,what is needed is a method and apparatus that efficiently deletes“large” boxes or parts of large boxes that using the interval Newtonmethod to solve the John conditions can only split.

SUMMARY

[0022] One embodiment of the present invention provides a system thatsolves a global optimization problem specified by a function | and a setof equality constraints q₁(x)=0 (i=1, . . . , r), wherein | is a scalarfunction of a vector x=(x₁, x₂, x₃, . . . x_(n)). During operation, thesystem receives a representation of the function | and the set ofequality constraints and stores the representation in a memory within acomputer system Next, the system applies term consistency to the set ofequality constraints over a subbox X, and excludes portions of thesubbox X that violate the set of equality constraints.

[0023] In a variation on the above embodiment, performing the intervalglobal optimization process involves analytically preconditioning theset of equality constraints through multiplication by an approximateinverse matrix B to produce a set of preconditioned equalityconstraints. It also involves applying term consistency to the set ofpreconditioned equality constraints over the subbox X, and excludingportions of the subbox X that violate any member of the set ofpreconditioned equality constraint,,.

[0024] In a variation on the above embodiment, performing the intervalglobal optimization process involves keeping track of a least upperbound |_bar of the function |(x), and removing from considerationsubboxes for which |(x)>|_bar. One way to do so involves applying termconsistency to the inequality |(x)<|_bar over the subbox X.

[0025] For any function of n-variables |(x) there are different ways toanalytically solve for a component x_(j) of the vector x. For example,one can write |(x)=g(x_(j))−h(x), where g(x_(j)) is a term in | forwhich it is possible to solve g(x_(j))=y for any x_(j) using g⁻¹(y). Foreach of these rearrangements, if a given interval box X is used as anargument of h, then the new interval X_(j) ⁺ for the j-th component ofX, is guaranteed to be at least as narrow as the original, X_(j).

X _(j) ⁺ =X _(j) ∩X′ _(j) where X′ _(j) =g ⁻¹(h(X)).

[0026] This process is then iterated using different terms g of thefunction |. After reducing any element X_(j) of the box X to X_(j) ⁺,the reduced value can be used in X thereafter to speed up the reductionprocess using other component functions and terms thereof if | is acomponent of the vector function |.

[0027] Hereafter, the notation g(x_(j)) for a term of the function |(x)implicitly represents any term of any component function. Thiseliminates the need for additional subscripts that do not add clarity tothe exposition.

[0028] In a variation on the above embodiment, applying term consistencyinvolves symbolically manipulating an equation within the computersystem to solve for a first term, g(x_(j)), thereby producing a modifiedequation g(x_(j))=h(x), wherein the first term g(x_(j)) can beanalytically inverted to produce an inverse function g⁻¹(y). Next, thesystem substitutes the subbox X into the modified equation to producethe equation g(X′_(j))=h(X), and solving for X′_(j) =g ⁻¹(h(X)).Finally, the system intersects X′_(j) with the interval X_(j) to producea new subbox X⁺, wherein the new subbox X⁺ contains all solutions of theequation within the subbox X, and wherein the size of the new subbox X⁺is less than or equal to the size of the subbox X.

[0029] In a variation on the above embodiment, performing the intervalglobal optimization process involves applying box consistency to the setof equality constraints q₁(x)=0 (i=1, . . . , r) over the subbox X, andexcluding portions of the subbox X that violate any member of the set ofequality constraints.

[0030] In a variation on the above embodiment, performing the intervalglobal optimization process involves evaluating a first terminationcondition. This first termination condition is TRUE if a function of thewidth of the subbox X is less than a pre-specified value, ε_(X), and theabsolute value of the function, |, over the subbox X is less than apre-specified value, ε_(F). If the first termination condition is TRUE,the system terminates further splitting of the subbox X.

BRIEF DESCRIPTION OF THE FIGURES

[0031]FIG. 1 illustrates a computer system in accordance with anembodiment of the present invention.

[0032]FIG. 2 illustrates the process of compiling and using code forinterval computations in accordance with an embodiment of the presentinvention.

[0033]FIG. 3 illustrates an arithmetic unit for interval computations inaccordance with an embodiment of the present invention.

[0034]FIG. 4 is a flow chart illustrating the process of performing aninterval computation in accordance with an embodiment of the presentinvention.

[0035]FIG. 5 illustrates four different interval operations inaccordance with an embodiment of the present invention.

[0036]FIG. 6 is a flow chart illustrating the process of finding aninterval solution to a nonlinear equation using term consistency inaccordance with an embodiment of the present invention.

[0037]FIG. 7A presents a first portion of a flow chart illustrating theprocess of using term consistency to solve an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention.

[0038]FIG. 7B presents a second portion of a flow chart illustrating theprocess of using term consistency to solve an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention.

[0039]FIG. 7C presents a third portion of a flow chart illustrating theprocess of using term consistency to solve an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention.

[0040]FIG. 7D presents a fourth portion of a flow chart illustrating theprocess of using term consistency to solve an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention.

[0041]FIG. 7E presents a fifth portion of a flow chart illustrating theprocess of using term consistency to solve an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention.

DETAILED DESCRIPTION

[0042] The following description is presented to enable any personskilled in the art to make and use the invention, and is provided in thecontext of a particular application and its requirements. Variousmodifications to the disclosed embodiments will be readily apparent tothose skilled in the art, and the general principles defined herein maybe applied to other embodiments and applications without departing fromthe spirit and scope of the present invention. Thus, the presentinvention is not limited to the embodiments shown, but is to be accordedthe widest scope consistent with the principles and features disclosedherein.

[0043] The data structures and code described in this detaileddescription are typically stored on a computer readable storage medium,which may be any device or medium that can store code and/or data foruse by a computer system. This includes, but is not limited to, magneticand optical storage devices such as disk drives, magnetic tape, CDs(compact discs) and DVDs (digital versatile discs or digital videodiscs), and computer instruction signals embodied in a transmissionmedium (with or without a carrier wave upon which the signals aremodulated). For example, the transmission medium may include acommunications network, such as the Internet.

[0044] Computer System

[0045]FIG. 1 illustrates a computer system 100 in accordance with anembodiment of the present invention. As illustrated in FIG. 1, computersystem 100 includes processor 102, which is coupled to a memory 112 anda to peripheral bus 110 through bridge 106. Bridge 106 can generallyinclude any type of circuitry for coupling components of computer system100 together.

[0046] Processor 102 can include any type of processor, including, butnot limited to, a microprocessor, a mainframe computer, a digital signalprocessor, a personal organizer, a device controller and a computationalengine within an appliance. Processor 102 includes an arithmetic unit104, which is capable of performing computational operations usingfloating-point numbers.

[0047] Processor 102 communicates with storage device 108 through bridge106 and peripheral bus 110. Storage device 108 can include any type ofnon-volatile storage device that can be coupled to a computer system.This includes, but is not limited to, magnetic, optical, andmagneto-optical storage devices, as well as storage devices based onflash memory and/or battery-backed up memory.

[0048] Processor 102 communicates with memory 112 through bridge 106.Memory 112 can include any type of memory that can store code and datafor execution by processor 102. As illustrated in FIG. 1, memory 112contains computational code for intervals 114. Computational code 114contains instructions for the interval operations to be performed onindividual operands, or interval values 115, which are also storedwithin memory 112. This computational code 114 and these interval values115 are described in more detail below with reference to FIGS. 2-5.

[0049] Note that although the present invention is described in thecontext of computer system 100 illustrated in FIG. 1, the presentinvention can generally operate on any type of computing device that canperform computations involving floating-point numbers. Hence, thepresent invention is not limited to the computer system 100 illustratedin FIG. 1.

[0050] Compiling and Using Interval Code

[0051]FIG. 2 illustrates the process of compiling and using code forinterval computations in accordance with an embodiment of the presentinvention. The system starts with source code 202, which specifies anumber of computational operations involving intervals. Source code 202passes through compiler 204, which converts source code 202 intoexecutable code form 206 for interval computations. Processor 102retrieves executable code 206 and uses it to control the operation ofarithmetic unit 104.

[0052] Processor 102 also retrieves interval values 115 from memory 112and passes these interval values 115 through arithmetic unit 104 toproduce results 212. Results 212 can also include interval values.

[0053] Note that the term “compilation” as used in this specification isto be construed broadly to include pre-compilation and just-in-timecompilation, as well as use of an interpreter that interpretsinstructions at run-time. Hence, the term “compiler” as used in thespecification and the claims refers to pre-compilers, just-in-timecompilers and interpreters.

[0054] Arithmetic Unit for Intervals

[0055]FIG. 3 illustrates arithmetic unit 104 for interval computationsin more detail accordance with an embodiment of the present invention.Details regarding the construction of such an arithmetic unit are wellknown in the art. For example, see U.S. Pat. Nos. 5,687,106 and6,044,454. Arithmetic unit 104 receives intervals 302 and 312 as inputsand produces interval 322 as an output.

[0056] In the embodiment illustrated in FIG. 3, interval 302 includes afirst floating-point number 304 representing a first endpoint ofinterval 302, and a second floating-point number 306 representing asecond endpoint of interval 302. Similarly, interval 312 includes afirst floating-point number 314 representing a first endpoint ofinterval 312, and a second floating-point number 316 representing asecond endpoint of interval 312. Also, the resulting interval 322includes a first floating-point number 324 representing a first endpointof interval 322, and a second floating-point number 326 representing asecond endpoint of interval 322.

[0057] Note that arithmetic unit 104 includes circuitry for performingthe interval operations that are outlined in FIG. 5. This circuitryenables the interval operations to be performed efficiently.

[0058] However, note that the present invention can also be applied tocomputing devices that do not include special.-purpose hardware forperforming interval operations. In such computing devices, compiler 204converts interval operations into a executable code that can be executedusing standard computational hardware that is not specially designed forinterval operations.

[0059]FIG. 4 is a flow chart illustrating the process of performing aninterval computation in accordance with an embodiment of the presentinvention. The system starts by receiving a representation of aninterval, such as first floating-point number 304 and secondfloating-point number 306 (step 402). Next, the system performs anarithmetic operation using the representation of the interval to producea result (step 404). The possibilities for this arithmetic operation aredescribed in more detail below with reference to FIG. 5.

[0060] Interval Operations

[0061]FIG. 5 illustrates four different interval operations inaccordance with an embodiment of the present invention. These intervaloperations operate on the intervals X and Y. The interval X includes twoendpoints,

[0062] x denotes the lower bound of X, and

[0063] x denotes the upper bound of X.

[0064] The interval X is a closed subset of the extended (including −∞and +∞) real numbers R* (see line 1 of FIG. 5). Similarly the interval Yalso has two endpoints and is a closed subset of the extended realnumbers R* (see line 2 of FIG. 5).

[0065] Note that an interval is a point or degenerate interval ifX=[x,x]. Also note that the left endpoint of an interior interval isalways less than or equal to the right endpoint. The set of extendedreal numbers, R* is the set of real numbers, R, extended with the twoideal points negative infinity and positive infinity:

R*=R∩{−∞}∪{+∞}.

[0066] In the equations that appear in FIG. 5, the up arrows and downarrows indicate the direction of rounding in the next and subsequentoperations. Directed rounding (up or down) is applied if the result of afloating-point operation is not machine-representable.

[0067] The addition operation X+Y adds the left endpoint of X to theleft endpoint of Y and rounds down to the nearest floating-point numberto produce a resulting left endpoint, and adds the right endpoint of Xto the right endpoint of Y and rounds up to the nearest floating-pointnumber to produce a resulting right endpoint.

[0068] Similarly, the subtraction operation X−Y subtracts the rightendpoint of Y from the left endpoint of X and rounds down to produce aresulting left endpoint, and subtracts the left endpoint of Y from theright endpoint of X and rounds up to produce a resulting right endpoint.

[0069] The multiplication operation selects the minimum value of fourdifferent terms (rounded down) to produce the resulting left endpoint.These terms are: the left endpoint of X multiplied by the left endpointof Y; the left endpoint of X multiplied by the right endpoint of Y; theright endpoint of X multiplied by the left endpoint of Y; and the rightendpoint of X multiplied by the right endpoint of Y. This multiplicationoperation additionally selects the maximum of the same four terms(rounded up) to produce the resulting right endpoint.

[0070] Similarly, the division operation selects the minimum of fourdifferent terms (rounded down) to produce the resulting left endpoint.These terms are: the left endpoint of X divided by the left endpoint ofY; the left endpoint of X divided by the right endpoint of Y; the rightendpoint of X divided by the left endpoint of Y; and the right endpointof X divided by the right endpoint of Y. This division operationadditionally selects the maximum of the same four terms (rounded up) toproduce the resulting right endpoint. For the special case where theinterval Y includes zero, X/Y is an exterior interval that isnevertheless contained in the interval R*.

[0071] Note that the result of any of these interval operations is theempty interval if either of the intervals, X or Y, are the emptyinterval. Also note, that in one embodiment of the present invention,extended interval operations never cause undefined outcomes, which arereferred to as “exceptions” in the IEEE 754 standard.

[0072] Term Consistency

[0073]FIG. 6 is a flow chart illustrating the process of solving anonlinear equation through interval arithmetic and term consistency inaccordance with an embodiment of the present invention. The systemstarts by receiving a representation of a nonlinear equation |(x)=0(step 602), as well as a representation of an initial box with X withX_(j) in X (step 604). Next, the system symbolically manipulates theequation |(x)=0 into a form g(x_(j))−h(x)=0, wherein the term g(x_(j))can be analytically inverted to produce an inverse function g⁻¹(y)(step606).

[0074] Next, the system substitutes the initial box X into h(X) toproduce the equation g(X′_(j))=h(X) (step 608). The system then solvesfor X′_(j)=g⁻¹(h(X)) (step 610). The resulting interval 7° C., is thenintersected with the initial interval X_(j) to produce a new intervalX_(j) ⁺ (step 612).

[0075] At this point, if X_(j) ⁺ is empty, the system can terminate.Otherwise, the system can perform further processing. This furtherprocessing involves setting temporarily saving X_(j) by settingX⁽⁰⁾=X_(j) and also, setting X_(j)=X_(j) ⁺ (step 614). Next, if w(X⁽⁰⁾)is sufficiently reduced at step 616, the system returns to step 606 foranother iteration of term consistency on another term g of |(x).Otherwise, the system terminates.

[0076] Examples of Applying Term Consistency

[0077] For example, suppose |(x)−x²−x+6. We can define g(x) x² andh(x)=x−6. Let X=[−10,10]. The procedural step is (X′)²=X−6=[−16,4].Since (X′)² must be non-negative, we replace this interval by [0,4].Solving for X′, we obtain X′=±[0,2]. In replacing the range of h(x)(i.e., [−16,4]) by non-negative values, we have excluded that part ofthe range h(x) that is not in the domain of g(x)=x².

[0078] Suppose that we reverse the roles of g and h and use theiterative step h(X′)=g(X). That is X′−6=X². We obtain X′=[6,106].Intersecting this result with the interval [−10,10], of interest, weobtain [6,10]. This interval excludes the set of values for which therange of g(X) is not in the intersection of the domain of h(X) with X.

[0079] Combining these results, we conclude that any solution of|(X)=g(X)−h(X)=0 that occurs in X=[−10,10] must be in both [−2,2] and[6,10]. Since these intervals are disjoint, there can be no solution in[−10,10].

[0080] In practice, if we already reduced the interval from [−10,10] to[−2,2] by solving for g, we use the narrower interval as input whensolving for h.

[0081] This example illustrates the fact that it can be advantageous tosolve a given equation for more than one of its terms. The order inwhich terms are chosen affects the efficiency. Unfortunately, it is notknown how best to choose the best order.

[0082] Also note that there can be many choices for g(x). For example,suppose we use term consistency to narrow the interval bound X on asolution of |(x)=ax⁴+bx+c=0. We can let g(x)=bx and computeX′=−(aX⁴+c)/b or we can let g(x)=ax⁴ and compute X′=±[−(bX+c)/a]^(¼). Wecan also separate x⁴ into x²*x² and solve for one of the factorsX′=±[−(bX+c)/(aX²)]^(½).

[0083] In the multidimensional case, we may solve for a term involvingmore than one variable. We then have a two-stage process. For example,suppose we solve for the term 1/(x+y) from the function|(xy)=1/(x+y)−h(x,y)=0. Let xεX=[1,2] and y ”Y=[0.5,2]. Suppose we findthat h(X,Y)=[0.5,1]. Then 1/(x+y) ε [0.5,1] so x+y ε [1,2]. Now wereplace y by Y=[0.5,2] and obtain the bound [−1,1.5] on X. Intersectingthis interval with the given bound X=[1,2] on x, we obtain the new boundX′=[1,1.5].

[0084] We can use X′ to get a new bound on h; but his may requireextensive computing if h is a complicated function; so suppose we donot. Suppose that we do, however, use this bound on oui intermediateresult x+y=[1,2]. Solving for y as [1,2]−X′, we obtain the bound[−0.5,1]. Intersecting this interval with Y, we obtain the new boundY′=[0.5,1] on y. Thus, we improve the bounds on both x and y by solvingfor a single term of |.

[0085] The point of these examples is to show that term consistency canbe used in many ways both alone and in combination with the intervalNewton algorithm to improve the efficiency with which roots of a singlenonlinear equation can be computed. The same is true for systems ofnonlinear equations.

[0086] Term Consistency and Equality Constrained Interval GlobalOptimization

[0087] FIGS. 7A-7E collectively present a flow chart illustrating theprocess of using term consistency in solving an interval globaloptimization problem with equality constraints in accordance with anembodiment of the present invention. Generally, we seek a solution in asingle box specified by the user. However, any number of boxes can byspecified.

[0088] The boxes can be disjoint or overlap. However, if they overlap, aminimum at a point that is common to more than one box is separatelyfound as a solution in each box containing it. In this case, computingeffort is wasted. If the user does not specify an initial box or boxes,we use a default box. The process described below finds the globalminimum in the set of points formed by the set of boxes. We assume theseinitial boxes are placed in a list L₁ of boxes to be processed.

[0089] Suppose the user of the process knows a point x_bar that isguaranteed to be feasible. If so, we use this point to compute aninitial upper bound |_bar on the global minimum |. If x_bar cannot berepresented exactly on the computer, the system forms a representableinterval vector X containing x_bar. We evaluate |(X) and obtain [lowerbound |(X), upper bound |(X)]. Even if rounding and/or dependence aresuch that X cannot be numerically proven to be certainly feasible, werely upon the user and assume that X contains a feasible point.Therefore, we set |_bar equal to the upper bound of |(X).

[0090] Also the user may know an upper bound |_bar on |* even though hemay not know where (or even if) | takes on such a value in the feasibleregion defined by the inequality constraints. If so, we set | bar equalto this known bound. If the known bound is not representable on thecomputer, the system rounds the value up to a larger value that isrepresentable.

[0091] If no feasible point is known and no upper bound on |* is known,we set |_bar=+∞. We assume the user has specified a box size toleranceε_(x) and a function width tolerance ε_(F).

[0092] The system specifies w_(R) and w₁. In doing so, the system setsw_(R)=0. If a single box X⁽⁰⁾ is given, the system sets w₁=w(X⁽⁰⁾),wherein w(X⁽⁰⁾) is the width of the box X⁽⁰⁾. If more than one box isgiven, the system sets w¹ equal to the width of the largest one.

[0093] Thus, to initialize our process, the system specifies ε_(X),ε_(F), w_(R), w₁, and the initial box(es). In addition, the systemspecifies a bound |_bar if one is known.

[0094] The steps of the process are to be performed in the order givenexcept as indicated by branching.

[0095] For each initial box X in the list L₁, the system evaluates |(X)(step 701).

[0096] Next, if |_bar<∞, the system deletes any box X from L₁ for whichthe lower bound of |(X)>|_bar. This can be done while applying termconsistency (step 702).

[0097] If L₁ is empty, the system goes to step 745. Otherwise, thesystem finds the box X in L₁ for which the lower bound of |(X) issmallest. For later reference, the system calls this box X⁽¹⁾. This boxis processed next. The system deletes X⁽¹⁾ from L₁ (step 703).

[0098] Next, if w[|/(X)]>ε_(F), the system goes to step 706 (step 704).

[0099] Then, if w(X)≦ε_(x), the system puts X in list L₂ and goes tostep 703 (step 705).

[0100] Next, the system applies term consistency to the constraintequations q₁,(x)=0 (i=1, . . . , r) for i=1, . . . , r. If it is provedthat there is no point in X that satisfies any one of the constraints,the system goes to step 703 (step 706).

[0101] If sup(|(X))≦|_bar, the system goes to step 711 (step 707).

[0102] Next, if w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in listL₂ and goes to step 703 (step 708).

[0103] Then, if X⁽¹⁾ (as defined in step 703) has been sufficientlyreduced, the system puts X in the list L₁ and goes to Step 703 (step709). We say that a box X is sufficiently reduced if any component of Xis reduced by an amount that is at least a fraction (say 0.25) of thewidth of the widest component of X.

[0104] If |[m(X)]≧|_bar, the system goes to step 713 (step 710).

[0105] For later reference, the system calls the current box X⁽²⁾. Thesystem then uses a linearization test to decide whether to linearize and“solve” the inequality |(x)≦|_bar. Let w₁ denote the width of thesmallest box for which M¹ is irregular. Let w_(R) denote the width ofthe largest box for which M¹ is regular. Initially, set w₁=w(X⁽⁰⁾) andw_(R)=0. We use w₁ and w_(R) to decide whether to generate and attemptto solve g(x)+J(x,X)(y−x)=0 for a given box X. We do so ifw(X)≦(w₁+w_(R))/2. (see Chapter 12 of Global Optimization Using IntervalAnalysis by Eldon R. Hansen, Marcel Dekker, Inc., 1992) (step 711).

[0106] If sup(|(X))≦|_bar, the system goes to step 714 (step 712).

[0107] Next, the system applies term consistency to the relation|(x)≦|_bar. If the result is empty, the system goes to step 703 (step713).

[0108] If w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂ andgoes to step 703 (step 714).

[0109] If X⁽¹⁾ (as defined in 703) has been sufficiently reduced, thesystem puts X in L₁ and goes to step 703 (step 715).

[0110] Next, the system applies box consistency to the constraintequations q₁(x)=0 (i=1, . . . , r) for i=1, . . . , r. If it is provedthat there is no point in X that satisfies any one of the constraints,the system goes go to step 703 (step 716).

[0111] If in |(|[m(X)])≧|_bar, the system goes to step 722 (step 717).

[0112] If the current box is the same box X⁽²⁾ defined in step 711, thesystem goes to step 720 (step 718).

[0113] Next, the system again uses the procedure in step 711 to try toreduce the upper bound |_bar (step 719).

[0114] If sup(|(X))≦f_bar, the system goes to step 722 (step 720).

[0115] Next, the system applies box consistency to the relation|(x)≦|_bar. If the result is empty, the system goes to step 703 (step721).

[0116] If w[(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂ andgoes to step 703 (step 722).

[0117] Next, if X⁽¹⁾ (as defined in step 703) has been sufficientlyreduced, the system puts X in the list L₁ and goes to step 703 (step723).

[0118] If in |(|[m(X)])≧f_bar, the system goes to step 713 (step 724).

[0119] Next, the system denotes the current box by X⁽³⁾, and applies alinear method involving a Taylor expansion to try to reduce X⁽³⁾ using|(X)≦|_bar. If the result is empty, the system goes to step 703 (step725).

[0120] Then, if w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in listL₂ and goes to step 703 (step 726).

[0121] If X⁽³⁾ as defined in step 726 was not sufficiently reduced inthe single step 725, the system goes to step 730.

[0122] Next, the system applies a quadratic method in which | isexpanded through quadratic terms to try to reduce the current box using|(X)≦|_bar. If the result is empty, the system goes to step 703 (step728).

[0123] If w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂ andgoes to step 703 (step 729).

[0124] If X⁽¹⁾ (as defined in step 703) has been sufficiently reduced,the system puts X in L₁ and goes to step 703 (step 730).

[0125] Next, the system uses a criterion w(X)>(w₁+w_(R))/2 to decidewhether to (try to) solve the constraints in linearized form. If not,the system goes to step 744 (step 731). Note that, we denotes the widthof the smallest box for which M¹=BJ(x,X) is irregular. If M¹ is regularfor a given box, w_(R) denotes the width of the largest box for which M¹is regular (step 731).

[0126] The system then replaces n-r of the variables by their intervalbounds, and renames the remaining variables as x₁l, . . . x_(r). Thenthe system linearizes the constraint functions as functions of thevariables now called x₁, . . . x_(r). The system also computes anapproximate inverse B of the approximate center of the Jacobian J(x,X)(step 732).

[0127] Next, the system preconditions the linearized system. If thepreconditioned coefficient matrix is regular, the system finds the hullof the linearized system. If the result is empty, the system goes tostep 703 (step 733).

[0128] If w[|/(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂and goes to step 703 (step 734).

[0129] Next, the system analytically multiplies the nonlinear system ofconstraint equations by the matrix B. It does so without replacing anyvariables by their interval bounds so that appropriate combinations andcancellations can be made. After the analytic multiplication iscomplete, the system replaces the fixed variables (as chosen in step732) by their interval bounds. The system also applies term consistencyto solve the i-th nonlinear equation of the preconditioned nonlinearsystem for the i-th variable (as renamed in step 732) for i=1, . . . ,r. If the result is empty, the system goes to step 703 (step 735).

[0130] If w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂ andgoes to step 703 (step 736).

[0131] Next, the system applies term consistency to solve the i-thnonlinear equation of the preconditioned nonlinear system for the i-th(renamed) variable for i=1, . . . , r. If the result is empty, thesystem goes to step 703. If the existence of a feasible point is proved(see Section 10.12), the system uses the result to update |_bar (step737).

[0132] If w[|(X)]≦ε_(F) and w(X)≦ε_(X) , the system puts X in list L₂and goes to step 703 (step 738).

[0133] Next, the system applies box consistency to solve the i-thnonlinear equation of the preconditioned nonlinear system for the i-th(renamed) variable for i=1, . . . , r. If the result is empty, thesystem goes to step 703 (step 739).

[0134] If w[|(X)]≦ε_(F) and w(X)≦ε_(X), the system puts X in list L₂ andgoes to step 703 (step 740).

[0135] Next, the system uses the criterion w(X)>(w₁+w_(R))/2 to decidewhether to (try to) solve the John conditions using a Newton method. Ifnot, the system goes to step 744 (step 741).

[0136] The system then applies one step of the interval Newton methodfor solving the John conditions. If the result is empty, the system goesto step 703. If the existence of a solution of the John conditions isproved, the system updates |_bar (step 742).

[0137] If the box X⁽¹⁾ (as defined in step 703) has been sufficientlyreduced, the system puts X in L₁ and goes to step 703 (step 743).

[0138] Any previous step that used term consistency, a Newton step, or aGauss-Seidel step might have generated gaps in the interval componentsof X. The system merges any such gaps when possible, and splits the box.This may involve deleting some gaps. The system places the subboxes(generated by splitting) in the list L₁ and goes to step 703 (step 744).

[0139] If the list L₂ is empty, the system indicates that there is nofeasible point in X⁽⁰⁾. The system then goes to step 750 (step 745).

[0140] If |_bar<+∞ and there is only one box in L₂, the system goes tostep 750 (step 746).

[0141] For each box X in L₂, if sup(|(m(X)))<|_bar, the system tries toprove existence of a feasible point. The system uses the results toupdate |_bar (step 747).

[0142] Next, the system deletes any box X from L₂ for which the lowerbound of f(X)>|_bar (step 748).

[0143] The system denotes the remaining boxes by X⁽¹⁾, . . . , X^((s))where s is the number of boxes remaining in L₂. The system thendetermines$\underset{\_}{F} = {{\min\limits_{1 \leq i \leq s}{{\underset{\_}{f}\left( X^{(1)} \right)}\quad {and}\quad \overset{\_}{F}}} = {\max\limits_{1 \leq i \leq s}{{\overset{\_}{f}\left( X^{(1)} \right)}\quad {\left( {{step}\quad 749} \right).}}}}$

[0144] Finally, the system terminates (step 750).

[0145] If, at termination, the list L₂ is empty, then all of the initialbox X⁽⁰⁾ has been eliminated. This provides proof that the initial boxX⁽⁰⁾ does not contain a feasible point.

[0146] Assume that at least one box remains in the list L₂. What we haveproved in this case depends on the final value of |_bar. If |_bar<+∞,then we know that a feasible point exists in the initial box X⁽⁰⁾. If|_bar =∞, there may or may not be a feasible point in X⁽⁰⁾.

[0147] Consider the case |_bar<∞. No matter how poor the bound |_bar on|*, we know that a global solution exists in X⁽⁰⁾; and it is in one ofthe remaining boxes. Also, we know that,

F≦|*≦{overscore (F)}.

[0148] If only one box X remains, then it must contain the globalsolution. In this case,

|(X)≦|*≦min {{overscore (|)}(X), |_(—) bar} and {overscore(|)}(X)−|(X)≦ε_(F).

[0149] Therefore, |(x)−|*≦ε_(F) for every point x in the box. Also,

x* ₁ , −X_(i) ≦ε _(X) {overscore (X)} ₁ −x* ₁≦ε_(X) (i=1, . . . n).

[0150] If more than one box remains, it is possible that one contains alocal solution at which | is less than our upper bound |_bar. Also,there might be more than one global solution occurring in separateboxes. We know only that

F≦|*≦min{{overscore (F)}, | _(—) bar}

[0151] and that the global minimum point(s) are in the remaining boxes.

[0152] If the final value of |_bar is ∞ X⁽⁰⁾ is not entirely deleted,then X⁽⁰⁾ may or may not contain a feasible point. We do not know. It ishighly probable that a solution exists since, otherwise, we expect allof X⁽⁰⁾ to be deleted. However, we do know that if a feasible point doesexist in X⁽⁰⁾, then

F≦|*≦{overscore (F)}

[0153] and x* is somewhere in the remaining box(es). All local solutionsin the initial box are contained in the final solution box(es).

[0154] It is possible that every point in the initial box X⁽⁰⁾ isinfeasible. However, our process can delete all of X⁽⁰⁾ (and thus provethere is no solution) only if every point in X⁽⁰⁾ is proved to befeasible (i.e., is certainly infeasible). Even if every point in X⁽⁰⁾ iscertainly infeasible, our process may still not delete all of X⁽⁰⁾. Theprobability of it doing so in this case is greater when the box sizetolerance ε_(X) is smaller.

[0155] Thus, we can prove there is no feasible point in X⁽⁰⁾, but we donot guarantee doing so when this is, in fact, the case.

[0156] Regardless of what procedure is used to delete a part of X⁽⁰⁾, weknow that the deleted part cannot contain a solution.

[0157] The foregoing descriptions of embodiments of the presentinvention have been presented only for purposes of illustration anddescription. They are not intended to be exhaustive or to limit thepresent invention to the forms disclosed.

[0158] For example, global optimization problems can contain bothinequality and equality constraints. In this case, the present inventioncan be combined with a system that solves an interval inequalityconstrained global optimization problem to form an overall system.

[0159] Accordingly, many modifications and variations will be apparentto practitioners skilled in the art. Additionally, the above disclosureis not intended to limit the present invention. The scope of the presentinvention is defined by the appended claims.

What is claimed is:
 1. A method for using a computer system to solve aglobal optimization problem specified by a function | and a set ofequality constraints, the method comprising: receiving a representationof the functions and the set of equality constraints q₁(x)=0 (i=1, . . ., r) at the computer system, wherein | is a scalar function of a vectorx=(x₁, x₂, x₃, . . . x_(n)); storing the representation in a memorywithin the computer system; performing an interval global optimizationprocess to compute guaranteed bounds on a globally minimum value of thefunction |(x) subject to the set of equality constraints; whereinperforming the interval global optimization process involves, applyingterm consistency to the set of equality constraints over a subbox X, andexcluding portions of the subbox X that can be shown to violate any ofthe equality constraints.
 2. The method of claim 1, wherein performingthe interval global optimization process involves: preconditioning theset of equality constraints through multiplication by an approximateinverse matrix B to produce a set of preconditioned equalityconstraints; applying term consistency to the set of preconditionedequality constraints over the subbox X; and excluding portions of thesubbox X that can be shown to violate any of the preconditioned equalityconstraints.
 3. The method of claim 1, wherein performing the intervalglobal optimization process involves: keeping track of a least upperbound |_bar of the function |(x); unconditionally removing fromconsideration any subbox for which in|(|x))>|_bar; applying termconsistency to the inequality |(x)≦f_bar over the subbox X; andexcluding portions of the subbox X that violate the inequality.
 4. Themethod of claim 1, wherein applying term consistency involves:symbolically manipulating an equation within the computer system tosolve for a term, g(x_(j)), thereby producing a modified equationg(x_(j))=h(x), wherein the term g(x_(j)) can be analytically inverted toproduce an inverse function g⁻¹(y); substituting the subbox X into themodified equation to produce the equation g(X′_(j))=h(X); solving forX′_(j)=g⁻¹(h(X)); and intersecting X′_(j) with the interval X_(j) toproduce a new subbox X⁺; wherein the new subbox X⁺ contains allsolutions of the equation within the subbox X, and wherein the size ofthe new subbox X⁺ is less than or equal to the size of the subbox X. 5.The method of claim, 1, wherein performing the interval globaloptimization process involves: applying box consistency to the set ofequality constraints q₁(x)=0 (i=1, . . . , r) over the subbox X; andexcluding portions of the subbox X that violate the set of equalityconstraints.
 6. The method of claim 1, wherein performing the intervalglobal optimization process involves: evaluating a first terminationcondition; wherein the first termination condition is TRUE if a functionof the width of the subbox X is less than a pre-specified value, ε_(x),and the absolute value of the function, |, over the subbox X is lessthan a pre-specified value, ε_(F); and if the first terminationcondition is TRUE, terminating further splitting of the subbox X.
 7. Themethod of claim 1, wherein performing the interval global optimizationprocess involves performing an interval Newton step on the Johnconditions.
 8. A computer-readable storage medium storing instructionsthat when executed by a computer system cause the computer system toperform a method for using a computer system to solve a globaloptimization problem specified by a function | and a set of equalityconstraints, the method comprising: receiving a representation of thefunction | and the set of equality constraints q₁(x)=0 (i=1, . . . , r)at the computer system, wherein | is a scalar function of a vectorx=(x₁, x₂, x₃, . . . x_(n)); storing the representation in a memorywithin the computer system; performing an interval global optimizationprocess to compute guaranteed bounds on a globally minimum vale of thefunction |(x) subject to the set of equality constraints; whereinperforming the interval global optimization process involves, applyingterm consistency to the set of equality constraints over a subbox X, andexcluding portions of the subbox X that can be shown to violate any ofthe equality constraints
 9. The computer-readable storage medium ofclaim 8, wherein performing the interval global optimization processinvolves: preconditioning the set of equality constraints throughmultiplication by an approximate inverse matrix B to produce a set ofpreconditioned equality constraints; applying term consistency to theset of preconditioned equality constraints over the subbox X; andexcluding portions of the subbox X that can be shown to violate any ofthe preconditioned equality constraints.
 10. The computer-readablestorage medium of claim 8, wherein performing the interval globaloptimization process involves: keeping track of a least upper bound|_bar of the function |(x); unconditionally removing from considerationany subbox for which in|(|(x))<|_bar; applying term consistency to theinequality |(x)≦|_bar over the subbox X; and excluding portions of thesubbox X that violate the inequality.
 11. The computer-readable storagemedium of claim 8, wherein applying term consistency involves:symbolically manipulating an equation within the computer system tosolve for a term, g(x_(j)), thereby producing a modified equationg(x_(j))=h(x), wherein the term g(x_(j)) can be analytically inverted toproduce an inverse function g⁻¹(y); substituting the subbox X into themodified equation to produce the equation g(X′_(j))=h(X); solving forX′_(j)=g⁻¹(h(X)); and intersecting X′_(j) with the interval X_(j) toproduce a new subbox X⁺; wherein the new subbox X⁻ contains allsolutions of the equation within the subbox X, and wherein the size ofthe new subbox X⁺ is less than or equal to the size of the subbox X. 12.The computer-readable storage medium of claim 8, wherein performing theinterval global optimization process involves: applying box consistencyto the set of equality constraints q₁(x)=0 (i=1, . . . , r) over thesubbox X; and excluding portions of the subbox X that violate the set ofequality constraints.
 13. The computer-readable storage medium of claim8, wherein performing the interval global optimization process involves:evaluating a first termination condition; wherein the first terminationcondition is TRUE if a function of the width of the subbox X is lessthan a pre-specified value, ε_(x), and the absolute value of thefunction, | over the subbox X is less than a pre-specified value, ε_(F);and if the first termination condition is TRUE, terminating furthersplitting of the subbox X.
 14. The computer-readable storage medium ofclaim 8, wherein performing the interval global optimization processinvolves performing an interval Newton step on the John conditions. 15.An apparatus that solves a global optimization problem specified by afunction | and a set of equality constraints, the apparatus comprising:a receiving mechanism that is configured to receive a representation ofthe function | and the set of equality constraints q₁(x)=0 (i=1, . . . ,r), wherein | is a scalar function of a vector x=(x₁, x₂, x₃, . . .x_(n)); a memory for storing the representation; an optimizer that isconfigured to perform an interval global optimization process to computeguaranteed bounds on a globally minimum value of the function |(x)subject to the set of equality constraints; wherein the optimizer isconfigured to, apply term consistency to the set of equality constraintsover a subbox X, and to exclude portions of the subbox X that can beshown to violate any of the equality constraints
 16. The apparatus ofclaim 15, wherein the optimizer is configured to: precondition the setof equality constraints through multiplication by an approximate inversematrix B to produce a set of preconditioned equality constraints; applyterm consistency to the set of preconditioned equality constraints overthe subbox X; and to exclude portions of the subbox X that can be shownto violate any of the preconditioned equality constraints
 17. Theapparatus of claim 15, wherein the optimizer is configured to: keeptrack of a least upper bound |_bar of the function |(x); unconditionallyremove from consideration any subbox for which in|(|(x))>|_bar; applyterm consistency to the inequality |(x)≦|_bar over the subbox X; and toexclude portions of the subbox X that violate the inequality.
 18. Theapparatus of claim 15, wherein while applying term consistency, theoptimizer is configured to: symbolically manipulate an equation to solvefor a term, g(x_(j)), thereby producing a modified equationg(x_(j))=h(x), wherein the term g(x_(j)) can be analytically inverted toproduce an inverse function g⁻¹(y); substitute the subbox X into themodified equation to produce the equation g(X′_(j))=h(X); solve forX′_(j)=g⁻¹(h(X)); and to intersect X′_(j) with the interval X_(j) toproduce a new subbox X⁺; wherein the new subbox X⁺ contains allsolutions of the equation within the subbox X, and wherein the size ofthe new subbox X⁺ is less than or equal to the size of the subbox X. 19.The apparatus of claim 15, wherein the optimizer is configured to: applybox consistency to the set of equality constraints q₁(x)=0 (i=1, . . . ,r) over the subbox X; and to exclude portions of the subbox X thatviolate the set of equality constraints.
 20. The apparatus of claim 15,wherein the optimizer is configured to: evaluate a first terminationcondition; wherein the first termination condition is TRUE if a functionof the width of the subbox X is less than a pre-specified value, ex, andthe absolute value of the function, |, over the subbox X is less than apre-specified value, ε_(F); and to terminate further splitting of thesubbox X if the first termination condition is TRUE.
 21. The apparatusof claim 15, wherein the optimizer is configured to perform an intervalNewton step on the John conditions.